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^ ', Abstract 

We describe the code hfodd which solves the nuclear Skyrme-Hartree-Fock problem 
by using the deformed Cartesian harmonic oscillator basis. The user has a possibility of 
choosing among various symmetries of the nuclear HF problem for rotating or non-rotating 
nuclei; they vary from the non-axial parity-invariant nuclear shapes, through those also 
breaking the intrinsic parity, towards the least-restrictive case corresponding to only one 
symmetry plane. The code provides a solution for a complete superdeformed rotational 
band in an v4~150 nucleus within one CPU hour of the CRAY C-90 supercomputer or 
' within two-three CPU hours of a fast workstation. 
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No. of bits in a word: 64 

Has the code been vectorised? : Yes 

No. of lines in distributed program: 19 438 (of which 8 354 are comments and separators) 

Keywords: Hartree-Fock; Skyrme interaction; Self-consistent mean-field; Nuclear many- 
body problem; Super deformation; Quadrupole deformation; Octupole deformation; Pair- 
ing; Nuclear radii; Single-particle spectra; Nuclear rotation; High-spin states; Moments 
of inertia; Level crossings; Harmonic oscillator; Coulomb field; Point symmetries 

Nature of physical problem 

The nuclear mean-field and an analysis of its symmetries in realistic cases are the main 
ingredients of a description of nuclear states. Within the Local Density Approximation, 
or for a zero-range velocity-dependent Skyrme interaction, the nuclear mean-field is local 
and velocity dependent. This allows an effective and fast solution of the self-consistent 
Hartree-Fock equations even for heavy nuclei and for different configurations, deforma- 
tions, excitation energies, or angular momenta. 

Method of solution 

The program uses the Cartesian harmonic oscillator basis to expand single-particle wave 
functions of neutrons and protons interacting by means of the Skyrme effective interac- 
tion. The expansion coefficients are determined by the iterative diagonalization of the 
mean field Hamiltonians or Routhians which depend nonlinearly on the local neutron and 
proton densities. Suitable constraints are used to obtain states corresponding to a given 
configuration, deformation or angular momentum. 

Restrictions on the complexity of the problem 

The main restriction is the CPU time required for calculations of heavy deformed nuclei 
and for a given precision required. One symmetry plane is assumed. Pairing correlations 
are only included in the BCS limit and for the conserved time-reversal symmetry. 

Typical running time 

One Hartree-Fock iteration for the super deformed, rotating, parity conserving state of 
1 g|Dy 86 takes about nine seconds on the CRAY C-90 computer. Starting from the Woods- 
Saxon wave functions, about fifty iterations are required to obtain the energy converged 
within the precision of about 0.1 keV. In case when every value of the angular velocity 
is converged separately, the complete superdeformed band with precisely determined dy- 
namical moments J^> can be obtained within one hour of CPU on the CRAY C-90, or 
within two to three hours of CPU on the SG Power Challenge L or IBM RS/6000 com- 
puters. This time can be often reduced by a factor of three when a self-consistent solution 
for a given rotational frequency is used as a starting point for a neighboring rotational 
frequency. 

Unusual features of the program 
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The user must have an access to the NAGLIB subroutine f02axe or to the ESSL subrou- 
tine zhpev which diagonalize complex hermitian matrices, or provide another subroutine 
which can perform such a task. 

The code is written in single-precision for the use on a 64-bit processor. The compiler 
option -r8 or +autodblpad (or equivalent) has to be used to promote all real and complex 
single-precision floating-point items to double precision when the code is used on a 32-bit 
machine. 

LONG WRITE UP 



1 Introduction 

The method of solving the Skyrme-Hartree-Fock (HF) equations in the Cartesian har- 
monic oscillator (HO) basis was described in the previous publication M which is referred 
to as I. The present paper is a long write up of the code hfodd which implements this 
method and solves the three-dimensional HF equations. A one- dimensional Skyrme-HF 
code restricted to the spherical symmetry has been published by P.-G. Reinhard 0. Al- 
though several HF codes allowing triaxial deformations exist || ^, [5|, ||, the code hfodd 
is the first one of that kind which is made available in the public domain. Several earlier 
beta versions of the code have already been distributed and used by various groups. The 
present version (vl.60n) is published for the first time and replaces all previous versions. 

In Sec. H we present the numerical tests of the code, and in Sees. ||, §, and [| we describe 
its input, output, and source files, respectively. 

2 Numerical tests 

Accuracy of the solution of the HF equations with the wave functions expanded onto 
the Cartesian HO basis, Eq. (1-70), i depends on the three parameters hu> x , hto y , hto z 
defining the HO frequencies in three Cartesian directions, and on the number M of the 
HO states included in the basis. In the code hfodd we use the standard prescription 
J7]. § to chose the HO states included in the basis, namely, the M states with the lowest 
HO single-particle energies, 

£n x n y n z = hu x (n x + ~) + huj y (n y + ~) + hu z (n z + ~), (1) 

are selected among those which have n x <N , n y <N , and n z <N , where N is the fixed 
maximum number of HO quanta. It should be noted that in general both M and Nq 
have to be specified to define the basis. Only for large N Q , the basis is defined solely by 
M and does not depend on Nq. In this case, the grid of points (n x ,n y ,n z ) defining the 
states included in the basis forms a pyramid in three dimensions, with the inclined face 
delimited by the condition e nxn nz <const. On the other hand, only for small values of Nq 
the basis is defined solely by Nq and does not depend on the energy cut-off. In this case 
the corresponding grid of points n x n y n z forms a cube of the size N Q . In all intermediate 

3 Symbol (1-70) refers to Eq. (70) of I. 
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cases the shape of the basis corresponds to a pyramid with the corners cut off, or to a 
cube with the corners cut off. Usually N is chosen large enough so that all the states 
allowed by the energy cut-off are included in the basis. The shape of the basis is printed 
in form of three projections by the code hfodd on the output file, see the TEST RUN 
OUTPUT below. 

In what follows we discuss only the results obtained with axially deformed bases for 
u±=uj x =ujy, and present the results in terms of the standard parameters 

luq = (uj^Uz) 1 / 3 and q = lu±/lu z . (2) 

(The code accepts of course in general a triaxially deformed HO basis what offers the user 
some additional flexibility in approaching the actual physical problem). In Ref. || the 
dependence of results on the basis parameters and on the basis size has been tested for 
iVo<10, which at that time was at the limit of the numerical feasibility. In the present 
study we discuss similar tests up to iVo=20 (spherical case) and iV =26 (deformed case). 

From the tests based on low values of Nq, the authors of Ref. [|J concluded that 
the basis parameters ujq and q have to be optimized for a given basis size in such a 
way that the HF energy is minimized. This procedure is consistent with the variational 
character of the HF theory in which the best approximation to the energy is obtained by 
a minimization of the energy. Although the variational principle does not ensure that the 
best approximation to other observables is obtained for the wave functions corresponding 
to the minimum of energy, this turned out to be the case for low values of N . The 
necessity to optimize the basis parameters was the real bottleneck of the method, because 
it requires many HF calculations to be performed before one physical solution can been 
obtained. 

The picture described above does not apply to calculations performed with relatively 
large values of N . Of course, the best energy is still obtained with optimized values of 
the basis parameters, but the energy depends now very weakly on the basis parameters 
and the error generated by using a non-optimal basis is small. Moreover, the best results 
for other observables are not necessarily obtained with the optimal basis parameters. 
This can be well understood by noticing that the minimum obtained from a very weak 
dependence of the energy on the basis parameters may be illconditioned, and may occur 
at rather accidental value. Therefore, for large values of Nq, the basis optimization does 
not bring any substantial improvement of precision. Instead of the optimization, suitable 
"physical" values of the basis parameters can be determined based on simple geometrical 
considerations. This is in line with results of Ref. ||, where it was shown that a partial 
optimization followed by the liquid-drop-model estimates reproduces the optimal basis 
parameters even for small values of Nq. 

In the code hfodd the physical basis parameters are determined as follows (cf. also 
Ref. II). First we define ujq by the standard value || multiplied by some phenomenological 
factor / close to unity, 

huj Q = f x 41MeV/A 1/3 . (3) 

Based on the experience gained from the diagonalization of the Woods-Saxon Hamiltonian 
on the HO basis [fT0| , |ITJ, the recommended value is f—1.2. This particular scaling of the 
oscillator frequencies turns out to help in stabilizing the results of the diagonalization 
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with respect to adding the new basis states. Second, suppose we wish to perform the HF 
calculations for a given set of multipole deformations a^, i.e., for a nuclear shape defined 
by the surface S 



£ : R(0,<j>) = c(a) 



1+ E E axiX^O,. 

A=0 fj,=-\ 



(4) 



where c(a) is a function of such that the volume enclosed by the surface £ does not 
depend on a. Due to the assumed y-simplex symmetry, see Sec. 3.1 of I, all multipole 
deformations are real, and only those with /x>0 are used in the code hfodd because 
then a A _ M =(-l)^a AA1 . 

The lengths of principal axes of the volume enclosed in the surface S can be defined 
as R x —R(tt/2,0), Ry=R(ir /2,tt/2), and R z —R(0, 0). For complicated shapes they do 
not necessarily correspond to the principal axes of the moment of inertia or quadrupole 
moment. However, when the quadrupole deformation dominates, this will roughly be so 
if the quadrupole deformation is considered in its intrinsic frame of reference, i.e., if one 
uses «2i=0. The relative values of the HO frequencies are then defined by the condition 

^xRx — UyRy = ^zRzi (5) 

while the common proportionality constant is given by the condition {u) x uj y u) z ) 1 ^=uJo. 

The above prescription is very helpful if the values of physical deformations are (at 
least roughly) known for a given nucleus. This is true in particular for the superde- 
formed shapes which have very similar deformations in a given region of nuclei. For 
example, the superdeformed nuclei in the ^4~150 region can be very well described in the 
basis corresponding to the 152 Dy nucleus with deformations a 2 o=0.61 and a 40 =0.10 |T2|. 



This gives the physical basis parameters, hu x =hcu y =6.2A6 MeV and huj z =11.200 MeV, i.e., 
^o=9.219 MeV and g=1.793. 

We may now present results of tests performed for different basis parameters and com- 
pare them with those obtained for the physical basis parameters. All calculations below 
have been performed for the Skyrme parametrization SkM* with the coupling constants 
Cf=C^=Cf^ s =0 in the energy functional (1-12). The latter conditions correspond to sim- 
plifications of the Skyrme functional usually assumed in studies made within the spatial 
coordinates |13|, [14]], and allow a comparison of our results with those obtained by an 



independent method and an independent HF code. 

From now on we are going to use the two adjectives introduced above and associated 
with the basis parameters: optimal, i.e., adjusted to a given basis size so as to minimize 
the HF energy, and physical, i.e., corresponding to the actual information about the 
deformations of nuclei studied, and calculated through Eqs. (||), @, and fl5|) based on the 
experimental result or earlier Strutinsky-type calculations. 



2.1 Spherical nuclei 

We begin with the results obtained for a spherical nucleus 208 Pb for which the spherical 
basis is used, q=l, and the results are studied in function of Hujq. The physical value for 
this nucleus is fru;o=8.304MeV, and is represented in Figs. |T]and0 by the vertical dashed 
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lines. The horizontal dashed lines represent very precise solutions obtained by the simple 
one-dimensional spherical code in spatial coordinates, and are considered as the exact 
values. The solid lines represent the results calculated by the code hfodd for N from 
10 to 18, and the asterisk represents the calculation with iVo=20. For the spherical bases 
the number of the HO states is given by M=(N +l)(N +2)(N +3)/6, and hence in Figs. 
| and H] it increases from 286 to 1771. 

The energies of 208 Pb, shown in Fig. [I, converge rather slowly to the exact value of 
— 1635.956 MeV. For iVo=10 the error is about 6 MeV and only for iVo>16 becomes smaller 
than 1 MeV. For every value of N , the dot denotes the minimum of energy which defines 
the optimal value of Huq. For iVo=10 and 12 the optimal values are rather remote from the 
physical value and, moreover, two minima appear for A^ =12. For iVo=14 the minimum 
occurs at the physical value of Huq, but for larger bases, iVo=16 and 18, the minimum 
becomes poorly defined and deviates from the physical value of huo. 

Fig. shows the corresponding results of calculations for the matter rms radii of 
208 Pb. For every curve the dots represent the minima of the energies shown in Fig. [I]. 
One should note that the vertical scale of the Figure is very much expanded and shows a 
narrow region of radii within about 0.3% around the exact value of _R rms =5.5546fm. For 
iVo=10 the radii significantly depend on huj but the value of -R rms at the physical basis 
parameter is very close to the exact result, while the value at the optimal basis parameter 
is much worse. For iVo=12 the value at the physical Kuoq is again much better than those 
at any of the two optimal bases. For iVo=14 the physical and optimal bases coincide and 
give the exact value with the error of 0.05%, while for larger bases the dependence of -R rms 
on Hujq becomes weak, and the results for the optimal bases are only slightly better than 
those for the physical basis. Altogether, for any value of iVo>10 the physical basis gives 
the radii precise up to ±0.003 fm, or up to ±0.05%. 



2.2 Deformed nuclei 

In order to study the properties of the optimization of the deformation of the basis, we 
have performed a series of calculations for the non-rotating superdeformed state in 152 Dy 
(hu=0) with hu fixed at the physical value, and for several different values of q and 
M. These results can be compared to those obtained in spatial coordinates for the same 
nucleus and force, and with the three-dimensional cubic grid of points spaced by Ax=1.0, 
0.8, or 0.7fm. In this case smaller spacing gives more accurate results, and one obtains 
r§, respectively, the total energies £=-1234.611, -1230.769, or -1230.104MeV and the 



proton quadrupole moments Q p =18.316, 18.369, or 18.376 b. 

These values of energies were obtained by using the finite difference expressions for 



derivatives, which are less precise than the Fourier expressions ||16|| . If the HF equations are 
solved using the first type of approximation but nevertheless at the end of the convergence 
the energies are calculated with the second method, one obtains a good and cost-effective 
estimate of the energy |16, 17j. Such a procedure yields |15] for Ax=0.7fm the value of 



£=—1229.365 which we may consider for the purpose of the discussion below as the exact 
result. This energy is shown by the horizontal dashed line in Fig. |3|. In practice, the 
three-dimensional calculations in the spatial coordinates are usually done with Ax=lfm 
and employ the finite difference method. For these standard conditions the binding energy 
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of 152 Dy is overestimated by about 5 MeV. 

Fig. ^| shows the total energies in 152 Dy calculated by the code hfodd for M=300, 
600, 900, or 1200. When q varies from 1.3 to 2.3 the maximum numbers of the HO 
quanta in the directions perpendicular/parallel to the elongation axis vary for M=300 
from 9/13 to 8/18 quanta. Similarly, for M=1200 and 9= 1.793 the HO basis contains 
14/26 quanta. Here and in the following, the results for M=300 are presented as broken 
lines with visible jumps corresponding to individual states entering and leaving the HO 
basis. Indeed, keeping a fixed number M of states with varying q means that a given 
HO orbital may cross the boundary e nx n y n z <const. For relatively small values of M this 
may create some effects on various quantities (visible in the expanded scales of presented 
Figures), while for larger values of M such effects disappear. 

The convergence of energy in function of M is again rather slow; for M=300 the error 
is around 5 MeV and it decreases to below lMeV only at M>600. For M=300 and 600 
the optimal basis deformations are rather different from the physical value. However, the 
energy gain from optimizing the value of q is very small. For M=900 the optimal and 
physical values of q are close to one another, but the minimum of energy is hardly visible. 
For M=1200 we obtain £—— 1229.383 in a very good agreement with the exact result. I 
should be noted that both estimates compared here are variational, and hence both may 
still slightly decrease with further improvement of the numerical parameters. 

Fig. f| shows the values of proton quadrupole moments calculated by the code hfodd 
and compared with the value of 18.376b obtained for Ax=0.7fm (see above). The latter 
value is for the purpose of the present test considered as the exact result. Again the region 
shown on the ordinate is very narrow; within 0.5% of the value of Q p . On the other hand, 
the region of q presented on the abscissa is fairly large; it corresponds to the standard 
quadrupole deformations f3=a.2o of the nuclear surface (J|) ranging from 0.3 to 1.0. 

For M=300 one obtains much a better value of Q p for the physical basis than for the 
optimal basis. This situation is reversed for M=600. For M=900 the dependence on 
basis parameters is very weak and the results are very close to the physical-basis value of 
<5 P =18.393 obtained for M=1200. The values of Q p calculated by the code hfodd tend 
to be slightly larger than the exact value obtained within the spatial-coordinate method. 
We have checked that the optimization of Tiloq does not remove this residual difference. 

Similarly as for tests performed in spherical nuclei, one can conclude that the optimiza- 
tion of the basis does not, in general, improve the results as compared to those obtained 
for the physical basis. Within a very limited space of M=300 HO states one may obtain 
a very reasonable estimate of Q p when the physical basis is used. Taking into account 
the fact that the discrepancy is larger for M=600, we may attribute a conservative error 
of 0.05 b, or 0.3% to the calculated values of the quadrupole moments. A precision of 
this order is perfectly sufficient for all practical purposes, and can of course be improved 
by using larger bases, instead of optimizing the basis parameters. In fact, the obtained 
precision is comparable to that of the spatial-coordinate method when it is used with the 
standard value of Ax=lfm (see above). 
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2.3 Rotating nuclei 

In Fig. H we show the values of the angular momentum / calculated for 152 Dy at hu— 
0.5 MeV. For all values of q and M considered in the Figure one obtains results which are 
within 0.1 h of the value calculated for M=1200. Again, there is no advantage in using 
the optimal bases as compared to the physical one. The optimal values of q were obtained 
by minimizing the Routhian (1-21) at constant hu. We do not show here the plot of the 
Routhian in function of q because such a plot is very similar to that of the energy at zero 
spin, Fig. |, except for an almost constant shift by — ul. 

The properties of superdeformed bands are very often described in terms of moments 
of inertia, which can be calculated as the static or dynamic (first or second) moments 

J (1) = -, (6a) 
u 

J<* = £, (6b) 

respectively. The precision of determining is therefore governed by that with which 
the total angular momentum I is calculated, and for /k<j=0.5MeV it is equal to about 
0.2h 2 /MeV. In practical calculations, the dynamic moment has to be calculated by 
the finite-difference expression 

J{ 2) = [M - ?M / 7 x 

u 2 — Ui 



which is precise up to 



SJ^-±(u 2 -u^ 



d 2 J^ 



24 y z 17 du 2 ' 

This expression for bj^ corresponds to the lowest-order term neglected in the Taylor 
expansion leading to the finite-difference formula (0). The curvature of J^ 2 \ which de- 
termines its precision 5J {2 \ is usually fairly small (200-300 fr 4 /MeV 3 at most), except 
at the crossing frequencies where the precise values of are physically not important 
anyhow. Therefore, the finite difference expression (^) can usually be used already at 
u 2 -u 1 ^0.1MeV/h. 

In Fig. H we present results for calculated at hu=0.55 MeV by using the values 
of / obtained at cc;i=0.5 and u 2 =0.Q MeV/h. In view of the precision of the total angular 
momentum /, estimated above to be equal to 0.1 h, one could have expected that the 
values of obtained in this way will be precise only up to about lh 2 /MeV. In fact, 
the real precision of J 1 ^ turns out to be an order of magnitude better. This is so because 
most of the errors related to the use of the finite HO basis cancel out, and the difference 
of angular momenta turns out to be much more precise that the total value of I itself. 
As a result, the dynamic moments shown in Fig. ^ are precise up to 0.1 h 2 /MeV. 
This observation concerns also many other observables, and will be discussed in the next 
section in more detail. 

2.4 Relative energies, moments, and alignments 

As already pointed out in Ref . || , the differences of energies (the relative energies) depend 
much less on the basis parameters than the total (absolute) energies. This is illustrated 
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in Fig. [7], where we show the proton separation energies S p in the SD 152 Dy nucleus, 
calculated at spin zero for different values of q and M. The values of S p have been 
calculated by subtracting the total energy of 152 Dy from the total energy of 151 Tb after 
both have been determined self-consistently by using the same HO basis. The SD state in 
151 Tb has been constructed by creating a hole in the 7r[651]3/2(+z) orbital, cf. Ref. |I2 |. 
As seen in the Figure, the separation energies are precise up to about 50keV, which is 
two orders of magnitude better than the total energies in either of 151 Tb or 152 Dy nucleus. 
Similar precision is also expected for the excitation energies in a given nucleus, i.e., for 
the differences of energies corresponding to different configurations or angular momenta. 

The same principle also works for the quadrupole polarizations [pl| which are obtained 
by subtracting the quadrupole moments of the SD states in 151 Tb and 152 Dy nuclei. This 
is shown in Fig. [8], where the obtained quadrupole polarizations are precise up to 0.01b, 
compared to the 0.05 b precision obtained for the values of proton quadrupole moments 
Q p discussed above. Similarly, the relative alignments between these two nuclei, i.e., 
the differences of angular momenta at a given angular frequency, hu=0.5MeV here, are 
presented in Fig. [| These values characterize the quality of the identical bands |T9|] , and 
have been discussed in Ref. [fT2 ]. Here we see that using the truncated HO basis they are 
precise up to 0.03 h, which is a factor three better than the values of the angular momenta 
discussed in Sec. |273. 



2.5 Summary 



Table [l] summarizes the results discussed in Sections |2j]-fZ4]. The column denoted "ab 



solute values" refers to the errors of observables calculated in 208 Pb or 152 Dy for M=300 
and for the physical basis parameters. The errors are deduced by comparing these results 
with the exact ones, where available, or with the results obtained with M=1200. The 
percentage errors are also given. The column denoted "relative values" refers to the errors 
of differences of observables calculated in 151 Tb and 152 Dy. Due to the fact that most 
of the errors are systematic (they originate from the truncated HO basis), they tend to 
cancel out, and the errors of relative observables are smaller than those of the absolute 
observables. Since the relative values are much smaller than the absolute values, the per- 
centage errors are larger in this case, but still only of the order of several per cent. In all 
cases the precisions of relative values are sufficient for a meaningful comparison with data 
or with other theories. Whenever more precise results are required, they can be easily 
obtained by increasing the size of the basis. The basis optimization is not a competing 
procedure which would lead to a better precision. 



3 Input data file 

Different aspects of the organization of the code hfodd are presented in this section 
together with the description of the input data file. The code reads the input data from 
the standard FORTRAN input file, prints results on the standard FORTRAN output file, 
and also writes two auxiliary files which below are called the review file and the record file. 
If required, the code may also start the iteration from the Woods-Saxon wave functions, 
which have to be provided in the specialized file called the Woods-Saxon file, or from the 
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previously recorded HF potentials provided in the file called the replay file. 
The general structure of the input data file is defined by the following rules: 

1. Input data file is an ASCII file composed of independent items. 

2. An item is composed of the keyword line, and of the data line which follows imme- 
diately the keyword line. Two items contain no data lines, see point ^| below. 

3. The keyword line contains in columns 1 through 10 a keyword which is a specific 
text defining the item. If the keyword has less than 10 characters it has to be 
obligatorily padded with appropriate number of trailing spaces. Columns beyond 
10 are ignored and can be used to place arbitrary comments or texts. 

4. The data line contains a prescribed number of data values which are read in the 
free FORTRAN FORMAT. The type of data (REAL, INTEGER, or CHARACTER) should 
match the definition of the given data line. The CHARACTER data should be provided 
starting from the 13-th column of the data line. 

5. Items can be separated by an arbitrary numbers of comment lines. A line is treated 
as a comment and ignored provided it does not contain in columns 1 through 10 
any of the valid keywords. 

6. Two items contain no data lines. The first one is defined by the keyword EXECUTE, 
and starts a calculation for the currently defined input parameters. The second one 
is defined by the keyword ALL_D0NE, and terminates the program. The latter item 
is not required if the compiler is able to properly recognize the end of input data 
file. When the end of input data file is encountered, the program terminates as if 
the keyword ALL_D0NE was found. 

7. A given item may occur any number of times in different places of the input data 
file. Only the last one of the same items occurring before a given EXECUTE keyword 
is taken into account. 

8. Between two consecutive items EXECUTE, the order of other items is arbitrary. These 
other items define the data set, i.e, the set of input parameters, for which the 
calculation will start at the moment when the next item EXECUTE is found. 

9. All input parameters have the default values defined in the code (subroutine pre- 
DEf). Therefore, if the input data file contains only the item EXECUTE, and no other 
lines, a calculation for the default values will be performed. In the present version, 
this will result in performing 50 iterations for the superdeformed state of 1 ggDy 86 at 
the angular frequency of hu=0.5 MeV. 

10. Only those items which define the values of input parameters which are different 
than the default values have to be included in the input data file. The values of 
input parameters defined by any item will stay in effect till another occurrence of 
the same item, or till the code terminates. This rule is valid irrespectively of how 
many EXECUTE items follow the given item. 
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Together with the FORTRAN source code in the file hf odd.f , two input data files are 
provided, dyl52-a.dat and dyl52-b.dat. The first one contains only one line with the 
keyword EXECUTE. The second one, which is reproduced in the section TEST RUN INPUT 
below, contains all valid keywords and the input parameters identical to the default values. 
Therefore, the results of executing the program for the both provided input data files are 
strictly identical. The second input data file serves only as a suitable pattern to modify 
the input parameters. However, it should be a good practice to include in the input data 
file only those items which modify the input parameters with respect to the default values. 
In this way, a danger of accidentally modifying some intricate numerical parameters of 
the code is minimized. 

In the following subsections we discuss all the valid keywords and their corresponding 
data lines. The default values are given in the data lines following the keyword lines. In 
the data lines we also show names of the corresponding variables which are referred to in 
the text. The meaning of several parameters that are introduced below in the context of 
their function in the code is explained in detail in I. 

3.1 General data 

Keyword: NUCLID 

86, 66 = IN_FIX, IZ_FIX 
Numbers of neutrons (IN_FIX) and protons (IZ_FIX) in the nucleus under consideration. 
Calculations for odd and odd-odd nuclei require IR0TAT=1. 

Keyword: ITERATIONS 

50 = NO ITER 

For N0ITER>0, the specified number of iterations is performed. The code starts with the 
iteration number by constructing the initial potentials. Specifying N0ITER=0 requests 
only this initial phase. If the iteration is restarted from the previously recorded potentials 
(IC0NTI=1), then the counting of iterations continues from the previous value on. 

3.2 Interaction 

Keyword: SKYRME-SET 

SKM* = SKYRME 

CHARACTER*4 acronym of the Skyrme force parameter set. Must start at the 13-th column 
of the data line. At present, valid acronyms are SKM*, SIII, SKP, and SKI1. Other sets 
of parameters can easily be included in the subroutine params. 

Keyword: EVE_SCA_TS 

1., 1., 1., 1., 1., 1, 1, 1., 1., 1., 1., 1. 

SRH0_T, SRH0_S,SRH0DT, SRHODS , SLPR_T, SLPR_S, 

STAU_T, STAtLS, 
SSCU_T, SSCILS, 
SDIV_T, SDIV_S 

By using this item in the data set the coupling constants corresponding to a given Skyrme 
parameter set can be arbitrarily scaled. This allows calculations with modified Skyrme 



11 



functionals [O]. The time-even coupling constants in the total-sum representation (I- 
14) are multiplied by the 12 numbers from SRHCLT to SDIV_S. The variables with names 
ending with _T and _S multiply the "total" and "sum" coupling constants, respectively. 
The variables with names containing the acronyms RHO, LPR, TAU, SCU, and DIV, multiply 
the coupling constants with superscripts p, Ap, r, J, and VJ, respectively, and those 
with RHOD multiply the density-dependent part of C p . Similar name convention is used 
for many other variables in the code hfodd. 

Keyword: 0DD_SCA_TS 

111111111111 
SSPI_T, SSPI_S,SSPIDT, SSPIDS,SLPS_T, SLPS_S, 

SCUR_T, SCUR_S, 
SKIS_T, SKIS_S, 
SR0T_T, SR0T_S 

Same as above but for the time-odd coupling constants. Acronyms SPI, LPS, CUR, KIS, and 
ROT correspond to coupling constants with superscripts s, As, T, j, and Vj, respectively, 
and those with SPID correspond to the density- dependent part of C s . 

Keyword: EVE_SCA_PM 

1., 1., 1, 1., 1., 1., 1., 1., 1., 1., 1., 1. 

SRHCLP, SRH0_M,SRH0DP, SRHODM , SLPR_P, SLPR_M, 

STAU_P, STAILM, 
SSCU_P, SSCU_M, 
SDIV_P, SDIV_M 

Same as above but for the time-even coupling constants in the isoscalar-isovector rep- 
resentation, Eq. (1-15). The variables with names ending with _P and _M multiply the 
isoscalar and isovector coupling constants, respectively. The total-sum scaling factors are 
used first, and the isoscalar-isovector scaling factors are used afterwards. 

Keyword: 0DD_SCA_PM 

1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1. 

SSPI_P, SSPI_M,SSPIDP, SSPIDM,SLPS_P, SLPS_M, 

SCUR_P, SCURJ4, 
SKIS_P, SKIS_M, 
SR0T_P, SR0T_M 

Same as above but for the time-odd coupling constants in the isoscalar-isovector repre- 
sentation, Eq. (1-15). 

Keyword: G -SCALING 

1.0, 1.0 = FACTGN, FACTGP 
For IPAIRI=1 the code hfodd solves the BCS equations with the neutron and proton 
pairing strengths defined in Ref. |2IJ. These values can be modified by defining here the 
multiplicative factors FACTGN and FACTGP for neutrons and protons, respectively. 

3.3 Symmetries 
Keyword: ROTATION 

1 = IROTAT 
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Calculation with time-reversal breaking will be performed for IR0TAT=1, while the time- 
reversal symmetry will be conserved for IR0TAT=0. In the latter case the calculations will 
be performed only for one value of the simplex, s—+i, which gives almost twice shorter 
execution times. IR0TAT=0 is incompatible with providing a non-zero value of the angular 
frequency or with attempting a calculation for an odd or odd-odd nucleus. 

Keyword: SIGNATUREY 

1 = ISIGNY 

Calculation with parity/signature symmetry conserved are performed for ISIGNY=1, while 
the parity (and signature) will be broken for ISIGNY=0. In the former case the parity 
blocks will be separately diagonalized and the densities will be summed only in the 1/8 
of all the Gauss-Hermite nodes. This may give almost twice shorter execution times. 
ISIGNY=1 requires ISIMTX=ISIMTZ=1, and ISIGNY=0 requires ISIMTX=0, or ISIMTZ=0, 
or ISIMTX=ISIMTZ=0, see Sec. 3.4 of I. 

Keyword: TSIMPLEXES 

1, 1 = ISIMTX, ISIMTZ 
Calculation with conserved symmetries given by a;-simplex T and/or ^-simplex 7 ", see Sec. 
3.4 of I, will be performed for ISIMTX=1 and/or ISIMTZ=1, respectively. These symme- 
tries will be broken for ISIMTX=0 and ISIMTZ=0. Values of ISIMTX and ISIMTZ should 
be compatible with that of ISIGNY in such a way that ISIGNY=ISIMTX*ISIMTZ. 

Keyword: PAIRING 

= IPAIRI 

Calculation without pairing correlations (pure HF) will be performed for IPAIRI=0 while 
the BCS pairing correlations will be included for IPAIRI=1. The latter case is incompat- 
ible with IR0TAT=1. 

3.4 Configurations 

Keyword: VACSIM_NEU 

43, 43 = KVASIM(O), KVASIM(l) 
Numbers of lowest neutron states occupied in the two blocks, denoted by (+) and (— ), 
of given simplexes, s—+i and s——i, respectively. These numbers define the simplex 
reference configuration from which the particle-hole excitations are counted. The defini- 
tions of simplex reference configuration and excitations are ignored unless ISIGNY=0 and 
IPAIRI=0. 

Keyword: VACSIM_PRO 

33, 33 = KVASIM(O), KVASIM(l) 

Same as above but for the numbers of proton states. 

Keyword: PHSIMP_NEU 

1, 00, 00, 00, 00 = NUPAH0, KPSIMP, KPSIMM, KHSIMP, KHSIMM 
Neutron particle-hole excitations in the simplex blocks. NUPAH0 is the consecutive number 
from 1 to 5 (up to five sets of excitations can be specified in separate items). Particles 
are removed from the KHSIMP-th state in the (+) block and from the KHSIMM-th state in 
the (-) block, and put in the KPSIMP-th state in the (+) block and in the KPSIMM-th 
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state in the (— ) block. At every stage of constructing excitations the Pauli exclusion 
principle has to be respected (particle removed from an occupied state and put in an 
empty state). Values equal zero have no effect. In practice, reasonable excitations can 
only be constructed by consulting the printed lists of single-particle states with their 
consecutive numbers in given blocks, see the TEST RUN OUTPUT below. 

Keyword: PHSIMP_PRO 

1, 00, 00, 00, 00 = NUPAH0, KPSIMP, KPSIMM, KHSIMP, KHSIMM 

Same as above but for the proton particle-hole excitations. 

Keyword: VACSIG_NEU 

22, 22, 21, 21 =KVASIG(0,0), KVASIG(0,1), 
KVASIG(1,0), KVASIG(1,1) 

Numbers of lowest neutron states occupied in the four blocks, denoted by (+,+), (+,—), 

(—,+), and (—,—), of given (paritysignature) combinations, i.e., (ir, r)—(+l,+i), 

(— and respectively. These numbers define the parity /signature reference 

configuration from which the particle-hole excitations are counted. The definitions of 

parity /signature reference configuration and excitations are ignored unless ISIGNY=1 and 

IPAIRI=0. 

Keyword: VACSIG_PR0 

16, 16, 17, 17 =KVASIG(0,0), KVASIG(O.l), 
KVASIG(1,0), KVASIG(l.l) 

Same as above but for the numbers of proton states. 

Keyword: PHSIGN_NEU 

1, 00,00,00,00, 00,00,00,00 =NUPAH0, 

KPPPSP, KPPPSM, KPPMSP, KPPMSM, 
KHPPSP, KHPPSM, KHPMSP, KHPMSM 
Neutron particle-hole excitations in the parity/signature blocks. Basic principles are the 
same as for the excitations in the simplex blocks. Particles are removed from the KHPPSP- 
th state in the (+,+) block, from the KHPPSM-th state in the (+,— ) block, from the 
KHPMSP-th state in the (-,+) block, and from the KHPMSM-th state in the (-,-) block, 
and put in the KPPPSP-th state in the (+,+) block, in the KPPPSM-th state in the (+,-) 
block, in the KPPMSP-th state in the (-,+) block, and in the KPPMSM-th state in the (-,-) 
block. 

Keyword: PHSIGN_PRO 

1, 00,00,00,00, 00,00,00,00 =NUPAH0, 

KPPPSP, KPPPSM, KPPMSP, KPPMSM, 
KHPPSP, KHPPSM, KHPMSP, KHPMSM 

Same as above but for the proton particle-hole excitations. 

3.5 Ensemble of specific parameters referred to as "Numerical data" 
Keyword: MAX_MULTIP 

2, 4, 4 = NMUC0N, NMUC0U, NMUPRI 

Maximum multipolarities A of multipole moments used in the code for the constraints, 
Eq. (1-22), surface term of the Coulomb field, Eq. (1-101), and printed on the output, 
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respectively. Values not larger than A=9 are currently allowed. In case of the conserved 
parity/signature, only even multipoles are used in the Coulomb field. 

Keyword: COULOMB 

80, 79, 0.25 = NUMCOU, NUMETA, FURMAX 
NUMCOU gives the number of points jV Coul , Eq. (1-100), used when summing up the Coulomb 
Green function. The dimensionless parameter d defining the size of the Coulomb paral- 
lelepiped, Eqs. (1-102) and (1-96), is given by NUMCOU* FURMAX. NUMETA gives the order of 
the Simpson integration of the solid harmonics on the faces of the parallelepiped, see Sec. 
5 of I. 

Keyword: SLOW-DOWN 

0.5, 0.5 = SL0WEV, SL0W0D 
The standard prescription to calculate the HF potential in the next iteration is to mix a 
given fraction e of the HF potentials from the previous iteration, with the fraction 1— e of 
potentials given by expressions (1-27). SL0WEV and SL0W0D give the values of e separately 
for the time-even and time-odd potentials. 

Keyword: EPS_HERMIT 

1.0E-14 = EPSHER 

Numerical precision requested for determining values of Hermite polynomials. 
Keyword: OPTI_GAUSS 

1 = I0PTGS 

For I0PTGS=1 and IREAWS=0, expression (1-94) is used to calculate the order of the 
Gauss-Hermite integration, and the input parameters NXHERM, NYHERM, and NZHERM are 
ignored. 

Keyword: GAUSHERMIT 

18, 18, 32 = NXHERM, NYHERM, NZHERM 
Orders L x , L y , and L z of the Gauss-Hermite integration in three Cartesian directions. 
Must be even. Ignored if I0PTGS=0 or IREAWS=1. 

3.6 Parameters of the HO basis 
Keyword: BASIS_SIZE 

15, 301, 800.0 = N0SCIL, NLIMIT, ENECUT 
The HO basis is composed of states having not more than iVo=NOSCIL quanta in either of 
the Cartesian directions, and not more than M=NLIMIT states in total, see Sec. ^ The 
states are added to the basis according to the increasing energy of the deformed harmonic 
oscillator (jl|). In case of degenerate HO states (e.g., for an axially deformed HO) the 
complete multiplets are included, so the actual number of states LDBASE can be slightly 
larger than NLIMIT. If NLIMIT>0, the cut-off energy ENECUT is ignored; otherwise all states 
having HO energy smaller then ENECUT are included in the basis, and NLIMIT is ignored. 

Keyword: H0MEGAZER0 

1.2 = FCH0M0 

The code uses the standard value of the spherical HO frequency huo=Al MeV / A 1 ^ 3 mul- 
tiplied by the scaling factor /=FCH0M0, see Eq. (0). 
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Keyword: SURFAC_PAR 

86, 66, 1.23 = INNUMB, IZNUMB, ROPARM 
The code hfodd calculates parameters of the HO basis, and the zero-iteration Nilsson 
potential, by defining the standard classical surface S, Eq. corresponding to the 
volume |tt*R0PARM 3 * ( I NNUMB+ IZNUMB) . 

Keyword: SURFACLDEF 

2, 0, 0.61 = LAMBDA, MIU, ALPHAR 
The code defines frequencies of the deformed HO in three directions by using relation (|^), 
uj x R x =u y Ry=uj z R z , where are the lengths of principal axes of the nuclear surface £ de- 
fined in the standard way (f|) by real deformation parameters a a^ALPHAR (LAMBDA, MIU). 
The overall factor is defined by (uj x tu y uj z ) 1 ^ 3 =uj . For details see Sec. ^. 

3.7 Constraints 
Keyword: OMEGAY 

0.5 = OMEGAY 

Value of the cranking rotational frequency c<jj=OMEGAY, Eqs. (1-23) and (1-24). Non-zero 
value requires IR0TAT=1. Pure linear constraint on the value of spin, uj y =uoj is used for 
IFLAGI=0. 

Keyword: MULTCONSTR 

2, 0, 0.01, 42.0, 1 = LAMBDA, MIU, STIFFQ, QASKED, IFLAGQ 
For IFLAGQ=1, the mass multipole moment of the given multipolarity A and fi is con- 
strained. Values of LAMBDA, MIU, STIFFQ, and QASKED correspond respectively to A, /i, 
C\ M , and in Eq. (1-22). For IFLAGQ=0, there is no constraint in the given multipo- 
larity. 

Keyword: SPINC0NSTR 

0. 0, 0.0, = STIFFI, ASKEDI, IFLAGI 

For IFLAGI=1, the quadratic constraint on spin is used together with the linear constraint. 
Values of STIFFI and ASKEDI correspond respectively to Cj and J y in Eq. (1-23). For 
IFLAGI=0, there is no quadratic constraint on spin, (but there still can be the standard 
cranking linear constraint defined by OMEGAY). 

3.8 Output-file parameters 
Keyword: PRINT- ITER 

1, 0, 1 = IPRSTA, IPRMID, IPRST0 

The code prints results for the first, middle, and/or last iteration if the corresponding 
parameters IPRSTA, IPRMID, and IPRST0 equal 1. 

Keyword: EALLMINMAX 

-12.0, 0.0 = EMINAL, EMAXAL 
The code prints tables of single-particle properties for states with values of the single- 
particle Routhians between EMINAL and EMAXAL. No table is printed unless EMINAL < 
EMAXAL. 
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3.9 Files 



Keyword: REVIEWFILE 

HFODD.REV = FILREV 
CHARACTER*68 file name of the review file. Must start at the 13-th column of the data 
line. The ASCII review file is written after calculating every data set (i.e., once per every 
EXECUTE item), provided IREVIE>0. The file is not rewound, so the results for several 
data sets can be accumulated in a single file. This is so provided the filename FILREV is 
not changed between the EXECUTE items. The file contains all relevant parameters and 
results of calculation in a form suitable for reading by another program. It is ment as an 
interface to programs which analyze and/or plot the results. The file contains sections 
defined by keywords (different than keywords used in the input data file described here). 
The detailed structure of the review file is not documented in the present write up, and 
can be inferred from inspecting the specimen produced by the sample run, and from the 
subroutine review. 

Keyword: REVIEW 

2 = IREVIE 

The review file will not be written if IREVIE=0. For IREVIE=2, the table of single-particle 
properties is included in the review file in addition to other results which are written for 
IREVIE=1. 

Keyword: RECORDFILE 

HFODD.REC = FILREC 
CHARACTER*68 file name of the record file. Must start at the 13-th column of the data 
line. The binary record file is written after each HF iteration. It contains complete 
information which allows restarting the iteration in another run of the code. To restart, 
one has to specify IC0NTI=1 and provide the name of the file by defining FILREP. In case 
of the computer crush, the record file contains the results of the last performed iteration. 
Upon a successful completion of the given input data set it contains the results of the 
last performed iteration. The file is always rewound before it is written, so the results for 
consecutive iterations do not pile up. 

Keyword: REPLAYFILE 

HFODD.REP = FILREP 
CHARACTER*68 file name of the replay file. Must start at the 13-th column of the data 
line. The binary replay file with the name defined in FILREP must exist if IC0NTI=1, 
and will be read. If the filenames FILREP and FILREC are identical, the replay file will be 
subsequently overwritten as a new record file. These feature is implemented to facilitate 
chaining of jobs which follow one another. 

Keyword: WOODSAFILE 

WOODS.WFN = FILWOO 
CHARACTER*68 file name of the Woods-Saxon file. Must start at the 13-th column of the 
data line. This file is read provided IREAWS=1. The binary Woods-Saxon file constitutes 
an interface between the Woods-Saxon code and the code hfodd. It contains the Woods- 
Saxon wave functions and numerous other parameters which define the current calculation. 
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The parameters read from the Woods-Saxon file overwrite the values provided in the input 
data file. Since the current version of the Woods-Saxon code will be published separately 
PHI , the feature of starting the iteration from the Woods-Saxon results is not documented 
in the present write up. 



3.10 Starting the iteration 



Keyword: RESTART 

= ICONTI 

For IC0NTI=1, results stored in the replay file (written as a record file in a previous 
run) will be used to start the iteration. The replay file name should be provided by 
defining FILREP. If the previous run was done with IREAWS=1, the current run must also 
use IREAWS=1, and the same Woods-Saxon file must be provided. This so because the 
Woods-Saxon file contains not only the information about the starting potential (which 
is ignored for IC0NTI=1) but also defines the HO basis. 

Keyword: READ_W00DS 

= IREAWS 

For IREAWS=1, the results stored by the previously performed Woods-Saxon calculation 
will be used to start the iteration and to define the HO basis. For IC0NTI=1 and IREAWS=1 
the Woods-Saxon file must also be provided, and will be used only to define the HO basis. 

Keyword: NILSSONPAR 

0, -1.175, -0.247, -1.175, -0.352, 11.17, 11.17, 6.28 = 

NILDAT, CNILSN, DNILSN, CNILSP, DNILSP, 

HBANIX, HBANIY, HBANIZ 
For IREAWS=0 and IC0NTI=0, the code starts the calculation from the Nilsson potential. 
If NILDAT=1, the Nilsson parameters C and D [Ref. ||22|| , Eq. (2.89)] are for neutrons given 
by CNILSN and DNILSN and for protons by CNILSP and DNILSP, while the HO deformation 
is defined by HBANIX, HBANIY, and HBANIZ. If NILDAT=0, the Nilsson parameters C and 
D are defined by Eq. (2.91) and Table 2.3 of Ref. ||22|| , while the HO deformation is the 
same as that of the HO basis defined in Sec. 3.6 . In the latter case the Nilsson parameters 
read from the input file are ignored. 



4 Output file 

Together with the FORTRAN source code in the file hfodd.f, an example of the output 
file is provided in dyl52.out. Selected lines from this file are presented in the section 
TEST RUN OUTPUT below. This output file corresponds to the input file dyl52-b . dat 
reproduced in the section TEST RUN INPUT below. Most of the information printed on 
the output file is self-explanatory. Here we only give some details which are not explicit 
in the output file. 

The output file begins with the information pertaining to the general parameters of 
the calculation, then gives information about the starting point of the iteration, provides 
the convergence report, and finally contains the results calculated at the last iteration. 
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Section CLASSICAL NUCLEAR SURFACE lists the deformation parameters used to define 
the nuclear surface, Eq. (|J), from which the basis parameters are derived, as described in 
Sec. |. It also gives the basis parameters such as the OSCILLATOR FREQUENCIES : HBAROX, 
HBAROY, and HBAROZ corresponding to hu x , hu y , and %uj z . 

Section PHYSICAL CONSTANTS gives the values of he in MeVfm (H_BARC), of hc/e 2 
(HBC0E2), of the neutron and proton masses in MeV/c 2 (XMASSN and XMASSP), of the 
kinetic-energy coefficient h 2 /2m before (HBMASS) and after the center-of-mass correction 
(HBMRPA), Eq. (1-10), of the elementary charge squared e 2 in MeVfm (ECHAR2), and of the 
coefficient preceding the integral in the Coulomb exchange energy (COULEX), Eq. (1-19). 
For details see the comments in the setbas subroutine. 

Section OSCILLATOR LENGTHS, CONSTANTS, and FREQUENCIES gives the values of 1/%, 
b^, and hu^, respectively, which characterize the HO basis in three Cartesian directions, 
H=x, y, z. 

Section BASIS CUT-OFF CONTROL PARAMETERS gives maximum numbers of the HO 
quanta in three directions NXMAXX, NYMAXX, and NZMAXX corresponding to N x , N y , and N z , 
Eq. (1-70), as well as the orders of the Gauss- Hermite quadratures NXHERM, NYHERM, and 
NZHERM, corresponding to L x , L y , and L z , Eq. (1-92). It also gives the number M of the 
HO states included in the basis as requested in the input file (NLIMIT) and as used in the 
calculation (LDBASE). 

Section SHAPE OF THE OSCILLATOR-BASIS DIAMOND gives the numbers of the HO 
quanta in a given direction for fixed numbers of the HO quanta in both remaining direc- 
tions. The output is arranged in such a way that the shape of the grid of points n x n y n z 
is clearly visualized by projections in every of the three directions. 

Section PARAMETER SET gives the name and the values of the Skyrme force parameters, 
while the following section COEFFICIENTS DEFINING THE SKYRME FUNCTIONAL gives the 
corresponding values of the coupling constants in the Skyrme functional (1-12). These 
values take into account the scaling factors, Sec. |3.2| , which are printed in the section 
SCALING FACTORS unless all are equal 1. 

Sections PARITY/SIGNATURE CONFIGURATIONS or SIMPLEX CONFIGURATIONS give the 
vacuum and particle-hole configurations requested in the input file for the case of the 
parity symmetry conserved or not conserved, respectively. 

Section CONVERGENCE REPORT gives the list of performed iterations. For each iteration 
one line is printed with the energy (1-9), stability (1-37), average values of the quadrupole 
moments Q 2 o an d Q22, total angular momentum, angular frequency (1-24), and the ratio 
of energies £/£, see Sec. 2.5 of I. The line corresponding to the last iteration is not printed 
because a more detailed information is given for the final state below this section. 

Section SINGLE-PARTICLE PROPERTIES lists the single-particle states calculated for 
the Nilsson, Woods-Saxon, or Hartree-Fock Routhian operators. For every state one 
line is printed which gives the value of the single-particle Routhian, cf. Eq. (31), the 
consecutive numbers in the parity/signature or simplex blocks, the quantum numbers 
[N,n z , A]Q of the asymptotic Nilsson state which has the largest component in the given 
state, the average value of the parity operator (in %), the average values of projections of 
the intrinsic and total angular momenta (in h), and their ratio called the g-factor. 

Sections MULTIPOLE MOMENTS give the values of neutron, proton, or mass multipole 
moments in units of (10fm) A . Traditional normalization factors are used as explained in 
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the comments in the subroutine DEFUNI. In particular, the A=0 moment corresponds to 
the number of particles. 

Sections ROOT-MEAN-SQUARE AND GEOMETRIC SIZES gives the rms average values of 
the radius and of the x, y, and z coordinates. In order to better visualize the size of the 



nucleus, the geometric sizes are also calculated by multiplying the rms radius by y5/3 

and the rms coordinates by y/E. 

Section DENSITY INTEGRALS IN THE SKYRME FUNCTIONAL gives the integrals of prod- 
ucts of densities, which appear in the Skyrme functional (1-12). The terms in the func- 
tional are identified by the acronyms described in Sec. |3.2j . 

Section CONTRIBUTIONS TO ENERGY IN THE SKYRME FUNCTIONAL gives the values of 
various terms which appear in the Skyrme functional (1-12). These contributions are the 
products of the coupling constants and of the density integrals described above. The sums 
of time-even and time-odd contributions are also printed. 

Section ANGULAR MOMENTA gives the average values of the total and intrinsic neutron, 
proton, and total angular momentum (in %). It also gives the corresponding values and 
contributions to the first moment of inertia J^—Iju (in h 2 /MeV). 

Sections NEUTRON CONFIGURATIONS and PROTON CONFIGURATIONS give a visual repre- 
sentation of states occupied in the parity /signature or simplex blocks. The lines denoted 
by CONF: give the configurations requested in the input data set, while those denoted 
by VACC: give the configurations characterizing the given HF state. By comparing the 
two sets one can verify whether the requested configuration has been obtained, and even- 
tually devise a new configuration to be calculated. The consecutive numbers printed 
in the horizontal direction correspond to the consecutive numbers in blocks printed in 
the section SINGLE-PARTICLE PROPERTIES. By comparing the two sections one can effec- 
tively associate the Nilsson labels with the calculated configurations and also prepare the 
configuration input data described in Sec. EO. 



Section ENERGIES gives a summary of the energies calculated for the HF state. The 
kinetic (1-10), single-particle (1-32), and pairing (1-20) energies are printed for neutrons 
(NEU), protons (PRO), and all particles (TOT). The direct (DIR) and exchange (EXC) 
Coulomb energies, Eqs. (I-17a) and (1-19), are printed together with their sum (TOT). 
The multipole (MULT) and cranking (SPIN) constraint energies, Eqs. (1-22) and (1-23), 
are printed together with the corresponding corrections (CORR. ) given by Eqs. (1-35) and 
(1-36), respectively. Then the rearrangement energy (1-34) is printed followed by the value 
of the Routhian (1-21) and the spin-orbit and Skyrme energies, the latter two split in the 
time-even (EVE) and time-odd (ODD) contributions. Finally, the total energies £ (1-33) 
and £ (1-9) are printed as (SP) and (FUN), respectively, while their difference (1-37) is 
printed as the stability (STAB). 

5 FORTRAN source file 

The FORTRAN source code in is provided in the file hfodd.f and can be modified in 
several places which are described in this section. 
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5.1 Dimensions of arrays 



The code HFODD uses the arrays' dimensions declared through the PARAMETER statements. 
This allows changing the dimensions and adapting the size of the reserved memory to the 
problem being solved. Whenever too small a dimension is defined the code aborts with 
a message indicating the dimension which should be increased. Substantial amount of 
memory is required only for arrays depending on the following PARAMETER values: 

PARAMETER (NDMAIN=16) 

Should be larger or equal to the input parameter NOSCIL defined in Sec. 
PARAMETER (NDBASE=307) 



Should be larger or equal to the input parameter NLIMIT defined in Sec. |3.6| . It should 
also be larger or equal to the actual size of the HO basis LDBASE, which can be larger 
than NLIMIT in case of degenerate HO states. 

PARAMETER (NDXHRM=19 , NDYHRM=19 , NDZHRM=33) 

Should be respectively larger or equal to the input parameters NXHERM, NYHERM, and 
NZHERM defined in Sec. gl| 

On vector machines, parameters NDXHRM, NDYHRM, NDZHRM, and NDBASE should be odd 
in order to minimize the risk of bank memory conflicts. By the same token, parameter 
NDMAIN should be even, because it defines the matrix dimensions beginning with 0. 



5.2 Vectorization properties 

As discussed in Sec 4.3 of I, the code hfodd has to operate by using seven-fold nested 
short loops, and this part does not perform well in a vector processor. However, it turns 
out that the loops can be artificially made longer in such a way that the final CPU time in 
a vector processor actually becomes much shorter. All the places where this trick has been 
applied can be identified in the source file by finding the lines beginning with CVECTOR 
and CSCALAR, for example: 

CVECTOR 

DO KZ=0 , LAZOXY (NX , NY) +LAZ0XY (MX , MY) 
CSCALAR DO KZ=NZ+MZ , , -2 

RESULT=RESULT+C0EF00 (KZ , NZ , MZ , 3) *T_AUXI (KZ) 
END DO 

The line beginning with CSCALAR should be made active on a scalar or superscalar 
machine. The line immediately below the line which begins with CVECTOR should be 
active on a vector machine. The results of calculation do not depend on which version of 
loops is activated. 



5.3 Library subroutines 

The code hfodd requires an external subroutine which diagonalizes complex hermitian 
matrices. In the present version, the code calls the NAGLIB subroutine f02axe. The 
call to this subroutine is performed through the interface subroutine zhpev, which can 
be replaced by the ESSL subroutine of the same name and the same calling parameters. 
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In the present version, the code also uses an external subroutine CGEMM from the 
BLAS CRAY library, which multiplies three complex matrices. By using the value of 
PARAMETER (I_CRAY=0) the call is diverted to the functionally equivalent subroutine 
zgemul, which is provided in the source file. However, the latter one can also be re- 
placed by the ESSL subroutine of the same name and the same calling parameters. 
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TEST RUN INPUT 



NUCLID 
ITERATIONS 

SKYRME-SET 
EVE_SCA_TS 
0DD_SCA_TS 
EVE_SCA_PM 
0DD_SCA_PM 
G.SCALING 

ROTATION 
SIGNATUREY 
TSIMPLEXES 
PAIRING 

VACSIM_NEU 
VACSIM_PRO 
PHSIMP_NEU 
PHSIMP_PRO 
VACSIG_NEU 
VACSIG_PRO 
PHSIGN_NEU 
PHSIGN_PRO 

MAX.MULTIP 

COULOMB 

SLOW-DOWN 

EPS.HERMIT 

OPTI_GAUSS 

GAUSHERMIT 



66 



50 



General data 



Interaction 



SKM* 

RHO 
1. 1. 

SPI 
1. 1. 

RHO 
1. 1. 

SPI 
1. 1. 

1.0 1.0 



RHOD LPR TAU SCU DIV 

1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 

SPID LPS CUR KIS ROT 

1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 

RHOD LPR TAU SCU DIV 

1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 

SPID LPS CUR KIS ROT 

1. 1. 1. 1. 1. 1. 1. 1. 1. 1. 



SIMP 

43 
SIMP 

33 



SIMM 

43 
SIMM 

33 



PART 
00 00 
PART 



2 

80 
0.5 



00 
PPSP 

22 
PPSP 

16 

00 
00 

4 
79 
0.5 



00 
PPSM 

22 
PPSM 

16 



PMSP 

21 
PMSP 

17 



PARTICLES 
00 00 
PARTICLES 
00 00 



Symmetries 



Configurations 



HOLE 
00 00 

HOLE 
00 00 
PMSM 

21 
PMSM 
17 



00 
00 



00 



00 

Numerical data 



HOLES 
00 00 

HOLES 
00 00 



00 
00 



1.00E-14 
1 



18 



18 



4 

0.25 



32 



BASIS.SIZE 










15 


301 


800. 


HOMEGAZERO 


1.2 






SURFAC_PAR 










86 


66 


1.23 


SURFAC_DEF 










2 





0.61 


SURFAC_DEF 










4 





0.10 



Parameters of the HO basis — 



OMEGAY 

MULTCONSTR 

SPINCONSTR 

PRINT-ITER 
EALLMINMAX 

REVIEWFILE 

REVIEW 

RECORDFILE 

REPLAYFILE 

WOODSAFILE 

RESTART 

READ_WOODS 

NILSSONPAR 

EXECUTE 
ALL.DONE 



Constraints 



0.50 
2 



-12. 



0.01 42. 1 

0.00 0. 

Output-file parameters 



Files 



HFODD . REV 

2 

HFODD. REC 
HFODD. REP 
WOODS. WFN 







-1.175 



Starting the iteration 

-0.247 -1.175 -0.352 11.170 11.170 
Calculate 

Terminate 



6.280 
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TEST RUN OUTPUT 



******************************************************************* 

* HFODD HFODD HFODD HFODD HFODD HFODD HFODD HFODD * 

* * 
******************************************************************************* 

* * 

* SKYRME-HARTREE-FOCK CODE VERSION 1.60N * 

* * 

* ONE SYMMETRY-PLANE AND NO TIME-REVERSAL SYMMETRY * 

* DEFORMED CARTESIAN HARMONIC-OSCILLATOR BASIS * 

* * 
******************************************************************************* 

* * 

* JACEK DOBACZEWSKI AND JERZY DUDEK * 

* * 

* CENTRE DE RECHERCHES NUCLEAIRES, STRASBOURG, 1993-96 * 

* * 
******************************************************************************* 

******************************************************************************* 

* CLASSICAL NUCLEAR SURFACE DEFINED FOR: N = 86 Z = 66 * 

* * 

* * 

* ALIO = ZERO ALII = ZERO * 

* AL20 = 0.610 AL21 = ZERO AL22 = ZERO * 

* * 

* AL30 = ZERO AL31 = ZERO AL32 = ZERO AL33 = ZERO * 

* * 

* AL40 = 0.100 AL41 = ZERO AL42 = ZERO AL43 = ZERO AL44 = ZERO * 

* * 
******************************************************************************* 

* HOMEGA= 9.2190 FCH0MO= 1.2000 * 

* * 

* OSCILLATOR FREQUENCIES: HBAROX= 11.1998 HBAROY= 11.1998 HBAROZ= 6.2464 * 

* * 

* MOMENTS OF INERTIA: XMOMFC= 90.2596 YMOMFC= 90.2596 ZMOMFC= 42.8300 * 

* * 

* CENTRES OF MASS: CMSXFC= 0.0000 CMSYFC= 0.0000 CMSZFC= 0.0000 * 

* * 
******************************************************************************* 

* PHYSICAL CONSTANTS: H_BARC=197 . 32891000 HBC0E2=137. 03602000 * 

* * 

* XMASSN=938. 90590000 XMASSP=938 . 27231000 * 

HBMASS= 20.73620941 HBMRPA= 20.59978698 * 

* ECHAR2= 1.43997841 COULEX= -1.06350868 * 

* * 
******************************************************************************* 

* * 

* OSCILLATOR LENGTHS: X= 1.9243099 Y= 1.9243099 Z= 2.5766958 * 

* OSCILLATOR CONSTANTS: X= 0.5196668 Y= 0.5196668 Z= 0.3880939 * 

* OSCILLATOR FREQUENCIES: X=ll. 1997759 Y=ll. 1997759 Z= 6.2464470 * 

* * 
******************************************************************************* 

* * 

* BASIS CUT-OFF CONTROL PARAMETERS: NXMAXX= 8 NYMAXX= 8 NZMAXX= 15 * 

* OPTIMUM NUMBERS OF GAUSS POINTS: NXHERM= 18 NYHERM= 18 NZHERM= 32 * 

* NLIMIT= 301 LDBASE= 306 MCOUNT=4096 * 

* ENECUT= 800.0000 ELIMIT= 112.7539 * 

* * 

t SHAPE OF THE OSCILLATOR-BASIS DIAMOND I 

* * 
******************************************************************************* 



* 




NZ 


===»> 





1 


2 


3 


4 


5 


6 


7 


8 


9 


10 


11 


12 


13 


14 


15 


* 




MAX. NX => 


8 


8 


7 


7 


6 


6 


5 


4 


4 


3 


3 


2 


2 


1 








* 
* 


NX= 





(15) I 


8 


8 


7 


7 


6 


6 


5 


4 


4 


3 


3 


2 


2 


1 








* 


NX= 


1 


(13) 1 


7 


7 


6 


6 


5 


5 


4 


3 


3 


2 


2 


1 


1 









* 


NX= 


2 


(12) I 


6 


6 


5 


5 


4 


4 


3 


2 


2 


1 


1 














* 


NX= 


3 


(10) 1 


5 


5 


4 


4 


3 


3 


2 


1 


1 


















* 


NX= 


4 


( 8) I 


4 


4 


3 


3 


2 


2 


1 






















* 


NX= 


5 


( 6) I 


3 


3 


2 


2 


1 


1 























* 


NX= 


6 


( 5) I 


2 


2 


1 


1 




























* 


NX= 


7 


( 3) I 


1 


1 
































* 


NX= 


8 


( 1) 1 





































^:fc)|E:+:>|E:+:*:+:>|c:+:)|t^>|E:4:^:4:>|E:+:)|E:+:^ 

* * 

* PARAMETER SET SKM* : T0= -2645.00 Tl= 410.00 T2= -135.00 T3= 15595.00 * 

* P0WER=0.1667 W=130 X0= 0.09000 Xl= 0.00000 X2= 0.00000 X3= 0.00000 * 

* - - * 
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******#*;+i*^***;+:*;+i*^***;+:#****^**^**^*^:**^^*sk**^**^***** 

* COEFFICIENTS DEFINING THE SKYRME FUNCTIONAL * 

* * 

******************************************************************************* 



* 






TOTAL (T) 


SUM(S) 


ISOSCALAR(P) 


ISOVECTOR(M) 


* 
* 


CRHO 




-1382 


012500 


780 


275000 


-991 


875000 


390.137500 


* 


CRHOD 




1299 


583333 


-649 


791667 


974 


687500 


-324.895833 


* 


CLPR 




-85 


312500 


34 


218750 


-68 


203125 


17.109375 


* 


CTAU 




68 


750000 


-68 


125000 


34 


687500 


-34.062500 


* 


CSCU 







000000 


68 


125000 


34 


062500 


34.062500 


* 


CDIV_ 




-65 


000000 


-65 


000000 


-97 


500000 


-32.500000 


* 
* 


CSPI 




-59 


512500 


661 


250000 


271 


112500 


330.625000 


* 


CSPID 







000000 


-649 


791667 


-324 


895833 


-324.895833 


* 


CLPS 







000000 


34 


218750 


17 


109375 


17.109375 


* 


CCUR 




-68 


750000 


68 


125000 


-34 


687500 


34.062500 


* 


CKIS 







000000 


-68 


125000 


-34 


062500 


-34.062500 


* 


CROT 




-65 


000000 


-65 


000000 


-97 


500000 


-32.500000 



******************************************** 
******************************************** 
* 

* PARITY/SIGNATURE CONFIGURATIONS : 
* 

* Y_*_5_!L!L?! 

* (++)~(+ : r( : +)~(— ) (++r~(+~-) 



I_C_L_E_S 

~~C-+T~C-~-) 



NEUTRONS : 
PROTONS : 



22 
16 



22 
16 



21 
17 



21 
17 



(++) 






H_0_L_E_S 

(^ : rc : +) 



(-) 






******************************************************************************* 

******************************************************************************* 



* CONVERGENCE REPORT 



ITER 



ENERGY 



-559.936804 

1 -1163.385944 

48 -1208.769024 

49 -1208.768731 



STABILITY 

-759.775974 
88.001030 
0.000468 
0.000400 



020 

54.178 
48.812 
41.807 
41.806 



022 

-0.027 
0.025 
0.067 
0.067 



SPIN OMEGA HOW NICE * 



98.777 
53.243 
49.581 
49.581 



0.500 0.424287 

0.500 1.081832 

0.500 1.000000 

0.500 1.000000 



** 
** 



************************************************************ 



* SINGLE-PARTICLE PROPERTIES: HARTREE-FOCK 



********** :11c****** 
***************** 



NEUTRONS 



******************************************************************************* 



* 


NO) 


ENERGY 


(++,+- 


-+ 


--) 


N 


nz 


A 


OMEG> 


<P> 




JY 




SY 


GFACT 


* 
* 


76) 


-11 


743 


( 0, 


21 


0) 


5 


3 


2 


3/2> 


-100 





103 


-0 


092 


-0.896 


* 


77) 


-11 


728 


( 0, 





20) 


5 


3 


2 


3/2> 


-100 





209 


-0 


101 


-0.485 


* 


78) 


-11 


433 


(19, 





0) 


4 


1, 


1 


l/2> 


100 





007 


-0 


182 


-3.E+01 


* 


79) 


-11 


334 


(0,19 





0) 


4 


1, 


3 


5/2> 


100 


-0 


129 


-0 


040 


0.307 


* 


80) 


-11 


241 


(20, 





0) 


4 


1, 


3 


5/2> 


100 


-0 


100 


-0 


189 


1.885 


* 


81) 


-11 


137 


(21, 





0) 


6 


5 


1 


l/2> 


100 


1 


119 


-0 


176 


-0.157 


* 


82) 


-11 


069 


( 0,20 





0) 


6 


5, 


1 


l/2> 


100 





984 





084 


0.086 


* 


83) 


-10 


881 


( 0,21 





0) 


4 


1, 


1 


l/2> 


100 


-0 


112 


-0 


201 


1.784 


* 


84) 


-10 


379 


(22, 





0) 


6 


4 


2 


5/2> 


100 


-0 


155 


-0 


064 


0.414 


* 


85) 


-10 


365 


( 0,22 





0) 


6 


4 


2 


5/2> 


100 





007 


-0 


024 


-3.300 


* 


86) 


-9 


557 


( 0, 





21) 


7 


6 


1 


3/2> 


-100 


2 


517 





041 


0.016 


* 


87) 


-7 


880 


( 0, 


22 


0) 


5 


2 


1 


3/2> 


-100 





903 





201 


0.222 


* 


88) 


-7 


840 


( 0, 





22) 


5 


2 


1 


3/2> 


-100 





442 





217 


0.490 


* 


89) 


-7 


738 


(23, 





0) 


4 





2 


5/2> 


100 


-0 


247 





161 


-0.650 


* 


90) 


-7 


732 


( 0,23 





0) 


4 





2 


5/2> 


100 


-0 


235 





162 


-0.690 


* 


91) 


-7 


562 


( 0, 


23 


0) 


5 


2 


1 


3/2> 


-100 


1 


068 





153 


0.144 


* 


92) 


-7 


221 


( 0, 





23) 


5 


1, 


4 


9/2> 


-100 


-0 


337 





028 


-0.083 


* 


93) 


-7 


221 


( 0, 


24 


0) 


5 


1, 


4 


9/2> 


-100 


-0 


337 





028 


-0.083 


* 


94) 


-7 


096 


( 0, 





24) 


7 


7 





l/2> 


-100 


-0 


134 


-0 


167 


1.252 


* 


95) 


-7 


076 


(0,24 





0) 


6 


4 





l/2> 


100 


1 


232 





328 


0.266 


* 


96) 


-6 


350 


( 0,25 





0) 


4 








l/2> 


100 





078 





455 


5.825 


* 


97) 


-6 


284 


(24, 





0) 


6 


3 


3 


7/2> 


100 





140 


-0 


044 


-0.313 



*** 
*** 

* DENSITY INTEGRALS IN THE SKYRME FUNCTIONAL * 

* * 
******************************************************************************* 



* 






TOTAL (T) 




SUM(S) 


ISOSCALAR(P) 


ISOVECTOR(M) 


* 
* 


DRHO 




17.566618 


8 


912767 


17.566618 


0.258917 


* 


DRHOD 




12.488132 


6 


334617 


12.488132 


0.181102 


* 


DLPR 




-3.719401 


-1 


886188 


-3.719401 


-0.052976 


* 


DTAU 




15.598603 


7 


983209 


15.598603 


0.367815 


* 


DSCU 




0.119929 





064074 


0.119929 


0.008219 


* 


DDIV_ 




0.823496 





419204 


0.823496 


0.014912 


* 
* 


DSPI 




0.018303 





010522 


0.018303 


0.002740 


* 


DSPID 




0.013147 





007544 


0.013147 


0.001941 


* 


DLPS 




-0.037372 


-0 


023155 


-0.037372 


-0.008937 


* 


DCUR 




0.060724 





031580 


0.060724 


0.002437 


* 


DKIS 




. 029460 





015277 


0.029460 


0.001094 


* 


DROT 




0.007844 





004322 


. 007844 


0.000799 
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* 
* 




CONTRIBUTIONS 


TO ENERGY IN THE 


SKYRME 


FUNCTIONAL 




* 
* 


******************************************************************************* 


* 




TOTAL (T) 


SUM(S) 


ISOSCALAR(P) 


ISOVECTOR(M) 


* 


* 


ERHO — 


-24277 


285107 


6954 


409417 


-17423 


888833 


101 


013143 


* 


* 


ERHOD = 


16229 


367796 


-4116 


181326 


12172 


025847 


-58 


839377 


* 


* 


ELPR = 


317 


311377 


-64 


543010 


253 


674755 


-0 


906388 


* 


* 


ETAU = 


1072 


403924 


-543 


856103 


541 


076525 


-12 


528704 


* 


* 


ESCU = 





000000 


4 


365047 


4 


085077 





279970 


* 


* 


EDIV_ = 


-53 


527218 


-27 


248236 


-80 


290827 


-0 


484627 


* 


* 
* 


SUM EVEN: 


-6711 


729227 


2206 


945789 


-4533 


317455 


28 


534017 


* 
* 


* 
* 


ESPI = 


-1 


089276 


6 


957479 


4 


962258 





905945 


* 
* 


* 


ESPID = 





000000 


-4 


901833 


-4 


271249 


-0 


630584 


* 


* 


ELPS = 





000000 


-0 


792320 


-0 


639413 


-0 


152907 


* 


* 


ECUR = 


-4 


174748 


2 


151393 


-2 


106350 





082995 


* 


* 


EKIS = 





000000 


-1 


040747 


-1 


003485 


-0 


037261 


* 


* 


EROT_ = 


-0 


509892 


-0 


280907 


-0 


764838 


-0 


025961 


* 


* 
* 
* 


SUM ODD: 


-5 


773917 


2 


093065 


-3 


823078 





142226 


* 
* 
* 



* MULTIPOLE MOMENTS IN UNITS OF (10 FERMI) **LAMBDA TOTAL * 

* * 
******************************************************************* 



* 


000 


=152 


0000 




* 
* 
* 


qio 




ZERO 


qn 


* 
* 


020 


= 41 


8059 


0,21 


* 


Q30 




ZERO 


031 


* 
* 


040 


= 4 


7902 


Q41 



ZERO 

ZERO Q22 = 0.0675 

ZERO Q32 = ZERO Q33 = 

ZERO Q42 = 0.0058 043 = 



ZERO 
ZERO 



Q44 



-0.0011 * 



******************************************************************************* 

* * 

* ROOT-MEAN-SQUARE AND GEOMETRIC SIZES IN FERMIS TOTAL * 

* R_RMS = 5.5420 X_RMS = 2.3844 Y_RMS = 2.3712 Z_RMS = 4.4052 * 

* * 

* R_GEO = 7.1547 X_GEO = 5.3318 Y_GEO = 5.3022 Z_GEO = 9.8504 * 

* * 
******************************************************************************* 
******************************************************************************* 

* * 

* ANGULAR MOMENTA AND THE FIRST MOMENTS OF INERTIA FOR OMEGA = 0.500000 MEV * 
******************************************************************************* 



NEUTRONS 
PROTONS 



ORBITAL 



28.21798 
19.29890 



SPINS 
INTRINSIC 



TOTAL 



1.21125 
0.85368 



29.42923 
20.15257 



ORBITAL 



56.43596 
38.59779 



J_Cl> 

INTRINSIC 

2.42250 
1.70736 



TOTAL * 



58.85846 
40.30515 



* TOTAL 47.51687 2.06493 49.58181 95.03375 4.12986 99.16361 * 

* * 
******************************************************************************* 

* NEUTRON CONFIGURATIONS * 



* 
* 




P 


S 


12 


13 


14 


15 


16 


17 


18 19 


20 21 


22 


23 


24 


25 


26 


27 


28 


29 


30 


31 


32 


* 
* 
* 


CONF: 


+ 


+ 


1 


1 


1 


1 


1 


1 


1 1 


1 1 


1 
































* 
* 


VACC: 


+ 


+ 


1 


1 


1 


1 


1 


1 


1 1 


1 1 


1 
































* 


CONF: 


+ 




1 


1 


1 


1 


1 


1 


1 1 


1 1 


1 
































* 
* 


VACC: 


+ 




1 


1 


1 


1 


1 


1 


1 1 


1 1 


1 
































* 


CONF: 




+ 


1 


1 


1 


1 


1 


1 


1 1 


1 1 



































* 
* 


VACC: 




+ 


1 


1 


1 


1 


1 


1 


1 1 


1 1 



































* 


CONF: 






1 


1 


1 


1 


1 


1 


1 1 


1 1 



































* 


VACC: 






1 


1 


1 


1 


1 


1 


1 1 


1 1 




































******************************************* 

* * 

* ENERGIES (MEV) * 

* * 
******************************************************************************* 



* 


KINETIC: 


(NEU)= 


1652 


730622 


(PRO)= 


1107 


885480 


(TOT)= 


2760 


616103 


* 


SUM EPS: 


(NEU)= 


-2034 


973323 


(PRO)= 


-1133 


181733 


(TOT)= 


-3168 


155056 


* 
* 


PAIRING : 


(NEU)= 





000000 


(PRO)= 





000000 


(TOT)= 





000000 


* 
* 


COULOMB : 


(DIR)= 


564 


434280 


(EXC)= 


-25 


354572 


(TOT)= 


539 


079708 


* 


CONSTR. 


(MULT)= 





000377 


SLOPE= 


-0 


003881 


CORR.= 


-0 


081133 


* 


CONSTR. 


(SPIN)= 


-24 


790903 


SLOPE= 





500000 


CORR.= 


-12 


395451 


* 
* 


REARRANGEMENT ENERGY FROM THE 


SKYRME DENSITY - 


-DEPENDENT 


TERMS= 


1009 


023720 


* 
* 


ROUTHIAN 


(TOTAL ENERGY PLUS MULTIPOLE 


AND SPIN CONSTRAINTS) = 


-1233 


559005 


* 


SPIN-ORB 


(EVE)= 


-80 


775454 


(ODD)= 


-0 


790799 


(TOT)= 


-81 


566253 


* 
* 


SKYRME: 


(EVE)= 


-4504 


783438 


(ODD)= 


-3 


680852 


(TOT)= 


-4508 


464290 


* 


TOTAL : 


(STAB)= 





000343 


(SP) = 


-1208 


768136 


(FUN)= 


-1208 


768479 



******************************************************************************* 
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Table 1: Precision of results obtained for various observables calculated by the code 
hfodd with the physical basis parameters and M=300. Absolute values correspond to 
the results for 152 Dy and 208 Pb, while the relative values correspond to the differences 
between 151 Tb and 152 Dy. The values denoted as "precision" represent also typical orders 
of magnitude of the accuracy obtained in other cases for the same numerical conditions. 





Absolute values 


Relative values 


Observable 


Precision % 


Precision % 


Energy £ 
Radius -R rms 
Quadrupole moment Q p 
Angular momentum / 
Static moment 
Dynamic moment 


5MeV 0.4% 
0.003 fm 0.05% 
0.05 b 0.3% 
0.1ft 0.2% 
0.2ft 2 /MeV 0.2% 
0.1h 2 /MeV 0.1% 


0.05 MeV 1% 

0.01b 1% 
0.03 ft 2% 
0.06ft 2 /MeV 2% 
0.1ft 2 /MeV 2-10% 
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Figure 1: Ground-state energies of 208 Pb (solid lines) calculated as functions of %ujq for 
several values of the number iVo of HO shells included in the basis. Dots denote the 
minima, the asterisk gives the value calculated at iVo=20 for the physical value of hu 
shown by the vertical dashed line. The horizontal dashed line represents the exact result. 



Figure 2: Same as in Fig. [T] but for the rms radii of 208 Pb. 



Figure 3: Energies of superdeformed state in 152 Dy (solid lines) calculated at hu=0 as 
functions of the basis deformation parameter q for several values of the number M of 
states included in the HO basis. Dots denote the minima, the asterisk gives the value 
calculated at M=1200 for the physical value of q shown by the vertical dashed line. The 
horizontal dashed line represents the exact result, see text. 



Figure 4: Same as in Fig. |] but for the proton quadrupole moments Q p of superdeformed 
152 Dy at hu=0. 



Figure 5: Angular momenta of the rotating superdeformed nucleus 152 Dy (solid lines) 
calculated at hu=0.5MeV as functions of the basis deformation parameter q for several 
values of the number M of states included in the HO basis. Dots denote the minima of 
the Routhian (1-21), the asterisk and the horizontal dashed line gives the value calculated 
at M=1200 for the physical value of q shown by the vertical dashed line. 



Figure 6: Same as in Fig. |5| but for the dynamic moment ([/J) of superdeformed 152 Dy 
at ftw=0.55MeV. 



Figure 7: Proton separation energies in SD 152 Dy (solid lines) relative to the SD state in 
151 Tb, i.e., the differences of energies £ SD ( 151 Tb)-£ SD ( 152 Dy), obtained at hu=0. Solid 
lines denote values calculated as functions of the basis deformation parameter q for several 
values of the number M of states included in the HO basis. Dots denote the minima of 
the energy in 152 Dy (Fig. |3|), the asterisk and the horizontal dashed line gives the value 
calculated at M=1200 for the physical value of q shown by the vertical dashed line. 
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Figure 8: Same as in Fig. |7| but for the quadrupole polarizations induced in the 152 Dy 
core by the 7r[651](3/2)(+i) hole, i.e., the differences of proton quadrupole moments 

Qf( 151 Tb)-Qf( 152 D V ). 



Figure 9: Same as in Fig. [7] but for the relative alignments between the 151 Tb and 152 Dy 
nuclei, i.e., the differences of angular momenta J( 151 Tb)— J( 152 Dy) at hu)=0.5 MeV. 
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